The Sun's Preferred Longitudes and the CoupUng of Magnetic 

Dynamo Modes. 
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ABSTRACT 

Observations show that solar activity is distributed non-axisymmetrically, concentrating at 
"preferred longitudes" . This indicates the important role of non-axisymmetric magnetic fields 
in the origin of solar activity. We investigate the generation of the non-axisymmetric fields and 
their coupling with axisymmetric solar magnetic field. Our kinematic generation (dynamo) model 
operating in a sphere includes solar differential rotation, which approximates the differential rota- 
tion obtained by inversion of helioseismic data, modelled distributions of the turbulent resistivity, 
non-axisymmetric mean helicity, and meridional circulation in the convection zone. We find that 
(1) the non-axisymmetric modes are localised near the base of the convection zone, where the for- 
mation of active regions starts, and at latitudes around 30°; (2) the coupling of non-axisymmetric 
and axisymmetric modes causes the non-axisymmetric mode to follow the solar cycle; the phase 
relations between the modes are found. (3) The rate of rotation of the first non-axisymmetric 
mode is close to that determined in the interplanetary space. 



Subject headings: Sun: magnetic fields — MHD 
1. Introduction 

Solar magnetic activity tends to cluster at "pre- 
ferred longitudes" (Vitinskij 1969 ; Bogart 1982; 
Bai 1987; Antonucci et al. 1990; Benevolenskaya et 
de Toma et al. 2000). Spacecraft data collected 
over three solar cycles indicate a persistent de- 
pendence of the solar wind speed and the radial 
component of the interplanetary magnetic field on 
solar longitude. Patterns are found to be rotat- 
ing with a fixed period (Neugebauer et al. 2000). 
This indicates the involvement of non-axisymmetric 
(i.e. longitude dependent) large-scale magnetic 
fields in the formation and evolution of solar activ- 
ity as seen in solar magnetograms (Ruzmaikin et al. 
Henney & Harvey 2003). 

The direct expansion of observed solar magne- 
tograms into spherical harmonics shows the pres- 
ence of axisymmetric (longitude-independent) and 
non-axisymmetric modes (Altschuler et al. 1974; 
de Toma et al. 2000). The axisymmetric mode 
dominates at solar minimum and the first non- 
axisymmetric mode is prevalent at solar max- 



imum (Ruzmaikin et al. 2001). Earlier studies 
(Stix 1974; Krause & Radler 1980; Ivanova & Ruzmaikin 1985) 
have shown that non-axisymmetric modes could 
be generated by the mean-field dynamo under 
al. 19!^j.^g^jjj constraints imposed on the distributions 
of the differential rotation and the mean kinetic 
helicity (a-effect). These earlier studies did not 
include the correct distribution, in radius and lat- 
itude, of the solar differential rotation curve, un- 
known at the time. A more recent model of Moss 
(1999), which includes the updated rotation curve, 
also displays the generation of a non-axisymmetric 
component of the field. 

The non-axisymmetric magnetic fields are gen- 
erated by axisymmetric sources but, in kine- 
matic mean-field dynamos, evolve independently 
of the axisymmetric field. In particular, the non- 
axisymmetric modes do not have the same period 
as the axisymmetric mode. However, the observa- 
tion that in the Sun the non-axisymmetric fields 
do have the same solar cycle as the axisymmetric 
fields, suggests a coupling between the modes that 
locks the oscillation periods of the two together. 
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Here we prove that this is in fact possible and 
it can be achieved by relaxing the assumption of 
axisymmetry for the a-effect. 

The coupling of different models by non-axisymmetric 
motions was first studied for the Earth's mag- 
netic fields by Bullard & Gellman (1954). Non- 
axisymmetry is a key ingredient in models of 
tachocline solar fltix tube dynamics that studied 
the instabilities leading to rise of the flux tubes 
(Choudhuri 1989; Fan ct al. 1993: Ferriz-Mas & Shiissler 
Caligari et al. 1995; 1998). In these models the 
formation of magnetic loops emerging at the so- 
lar surface as active regions is considered to be 
a consequence of these non-axisymmetric (kink- 
type) instabilities. Although the dynamo itself is 
not treated in these studies, longitudinal insta- 
bilities are connected to the emergence of active 
regions, see also Ruzmaikin (2001). A calcula- 
tion of the a-effect from these instabilities, have 
been performed by (Ferriz-Mas et al. 1994). Such 
an a-effect, however, has not yet been included 
in mean- field dynamo models. Kinetic helicity 
in hydrodynamic simulation of shallow water in- 
stabilities used in models of the tachocline, also 
reveal a non-axisymmetric, to — 1, pattern for the 
a-effect (Dikpati & Oilman 2001). That is not in- 
cluded, however, in the dynamo model that the 
same authors have analysed. A non-axisymmetric 
a-effect has been considered in mean-field dynamo 
models of the galactic dynamo (Moss et al. 1991). 

Here for the first time we study the solar mean- 
field dynamo with a non-axisymmetric a. We dis- 
cuss the distribution in radius latitude and az- 
imuth of the first non-axisymmetric mode and re- 
late it to preferred longitudes and clustering of 
magnetic activity. It is important to note that, 
contrary to the axisymmetric fields, the non ax- 
isymmetric modes are rotating structures in the 
frame of reference of the body. The rate of rota- 
tion of these modes is determined by the solution 
of the dynamo equations. We determine the rate 
for the TO = 1 mode by looking at the evolution 
of the surface patterns of its radial magnetic field. 
Because the surface poloidal field expands into the 
solar wind, we relate this rate to the signatures of 
preferred longitudes found in the solar wind. 

Our kinematic model includes an approxima- 
tion to the helioseismically observed distribution 
of the differential rotation and modelled distribu- 
tions of other sources of the mean field generation. 



2. The Non- Axisymmetric Dynamo Model 

When the motions arc given, the generation of 
the mean magnetic field is described by the equa- 
tions 

a^B = V X (u X B-l-aB -r;B), (1) 

V • B = (2) 

(Mofffat 1978; Parker 1979; Krause & Radler 1980), 
^^^e B is the mean magnetic field (the first sta- 
tistical moment of the total field) , u is the mean 
velocity, a is the mean-helicity, and rj is the elec- 
trical resistivity, which includes the turbulent re- 
sistivity. 

Solar rotation is approximated by a simple an- 
alytical fit (Charbonneau et al. 1999; Moss 1999) 
to the distribution reconstructed by hclioscismic 
methods (Kosovichev 1996; Howe et al. 2000), 
see Fig. 1. Details are given in the Appendix. 
The rate of the equatorial rotation is taken as 
neg/27r = 460.7 nHz = 4.607 x lO-^s"!. This dis- 
tribution matches the core rotation Qq smoothly 
in a layer of thickness 0.04i?Q at O.69R0. The 
core rotation is about equal to the surface rota- 
tion fig at 30° latitude. The surface shear layer 
is not included in this model. However, the deep 
shear layer at O.7R0, i.e. the tachocline, is well 
approximated. 

The distribution of a cannot yet be deter- 
mined from observations. We therefore con- 
sider three typical cases, see Appendix. The 
radial profiles are shown in Fig. 2. In the first 
case. Ml, the distribution of a is independent 
of radius in the convection zone and rapidly de- 
creases below the convection zone;. In the sec- 
ond case, M2, a is concentrated near the so- 
lar surface and does not overlap with the shear 
layer. In the third case, M3, a is localised above 
the tachocline (Oharbonneau & MacOregor 1997; 
Mason et al. 2002) overlapping with the radial 
shear. The latitudinal dependence of a. in all 
three cases is assumed to be cos^. This simple 
angular dependence is partially justified by the 
fact that the pseudo-scalar a can be formed from 
the scalar product of the density gradient and 
rotation: a ~ Vp • ~ cos^. More elaborate 
angular dependences are produced by different 
models of a, see for example Riidiger & Branden- 
burg (1995), Dikpati & Oilman (2001), and can 
easily be incorporated into our model. 
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The resistivity rj is taken constant throughout 
the convection zone with a sharp decrease in the 
transition layer to the highly conductive radiative 
interior. We consider the case when the transition 
layer is coincident with the rotational shear layer. 
For computational needs, we assume that the re- 
sistivity near the base of the convection zone is 
50 to 200 times less than the turbulent resistivity 
in the bulk of the convection zone. (A realistic 
number would be of the order of 10^). We use 
Rq/t] as the time unit. We assign the value of 
the resistivity r] by considering the marginal solu- 
tion of model Ml. We match the equatorial rota- 
tion rate in Fig. 1 to the observed one, obtaining 
?7o = 2.74 X lO^^cm^s"^. Estimates made using 
classical mixing length theory, give a value for the 
resistivity of order lO^^cm^s""'^. 

The ratio Tc/Trot of the cycle period to the 
rotation rate of the body is non-dimensional and 
does not therefore depend on the choice of r?. For 
the Sun, the ratio of the cycle and the equato- 
rial rotation rate, is approximately 150. Higher 
values of this number imply a cycle period of the 
solutions longer than the solar one. Shorter val- 
ues, a shorter cycle period. We characterise the 
cycle periods of our solution in terms of this non- 
dimensional ratio. 

In the case of non-axisymmetric dynamo an 
analogous ratio of the rate of rotation of the 
non-axisymmetric mode to the body rotation rate 
could also be used to characterise the solutions. 
Such a number clearly does not exist in the case of 
the axisymmetric dynamo where no rotating struc- 
tures exist. 

2.1. Mathematical Formulation and Nu- 
merical Procedure 

We have developed a new code that solves the 
kinematic dynamo equations (1,2) in a spherical 
domain. We define a regular grid in the meridional 
circular sector r G [rj„,i?0],^^ G [0, tt], where 9 is 
the co-latitude. In the azimuthal direction cj) we 
expand all functions in m-modes, e*™"^. Solving 
the equations in the real space has advantages over 
using the spectral Legendre Transform. The m- 
modal transform can be numerically performed in 
a very efficient way by means of the Fast Fourier 
Transform algorithm. No such algorithm exists for 
the Legendre Transform, needed for the latitudinal 
direction 9. This makes spherical spectral codes 



numerically expensive to run. It is also easier to 
parallelise a code in a regular domain. 

We represent the divergence-free magnetic field 
by two scalar functions T and S, called the toroidal 
and poloidal potentials 

B = -r X VT-KV X (-r X VS). (3) 
^ « ' ^ ^ ' 

Toroidal Poloidal 

(Mofffat 1978; Krause & Radler 1980). This rep- 
resentation ensures that the field is divergence 
free. The curl of a poloidal field is a toroidal field 
and vice versa. 

The choice of T and S is made unique by means 
of the gauge condition 

/ T sin6» d6l # = , / S sin 9 d9 dcp = 0. 

(4) 

When the field vanishes, B = 0, then its toroidal 
and poloidal components Bp and Bf vanish as 
well. The vector equation (1) therefore reduces to 
two coupled equations for the toroidal and poloidal 
potentials. 

When a = a{r,9) and r] = rj{r), the governing 
Eqs (1), in spherical coordinates r, 6^ (/), reduce to 

dtT = RnVn + RaVa + RmVm 

+ f]\7^T + drT]--dr(rT), 
r 

dtS = RnUn + Ralla + RmUm 

+ TiV^S (5) 

with the gauge conditions (4) . The non-dimensional 
numbers Ra, Rm are defined as: 

rio-R© _ UqRq _ UmRq 
Hn — , -fta — , JIM — • 

Vo Vo Vo 

(6) 

Subscripted quantities indicate typical ampli- 
tudes. The functions Uq, Vq, Ua, Va, Um, Vm, are 
related to the toroidal and poloidal parts of the 
sources: 

{{n X r) X B)t = -r X VC/n, (7) 
V X X r) X B)p = -rxVVh, (8) 
(uM X B)t = -r X VUm, (9) 
Vx(umxB)p = -rxVVM, (10) 
(aB)T = -rxVC/„, (11) 
Vx(q!B)p = -rxVVa- (12) 
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Wc have split the velocity field u into its rotational 
(toroidal) and meridional components 

u{r,0) = n{r,9)xr + UM{r,e). (13) 

The functions Uf,Ve,Ua,Va,UM,VM can be 
calculated from the above relations: 

-LVo = r • V X V X X r) X B), (14) 

^L'^Un = r • V X ((J7 X r) X B), (15) 

-LVm = r • V X V X (UM X B), (16) 

-L'^Um = r-Vx(uMxB), (17) 

-L'^Va = r • V X V X (a(r,6') B), (18) 

-L'^Uc = r-V x{a{r,e)B), (19) 

where 

sm 86 \ 86 J sm 6 o(jr 
is the angular part of the Laplacian operator 

J.2 Qj. Qj. J.2 

The inversions in Eqs. (14-19) can be efficiently 
carried out numerically. 

2.2. The m-Modal Expansion 

We expand all functions T, S, V, U in Eq. (5) 
into the m-modes 

N 

T{r, 9,(1)) T-^{r, e)e''^* + cc, ... (22) 

m=0 

where T"^{r,6) are complex functions and N is 
the maximum order of the expansion. When the 
sources (u, a, rj) are axisymmetric, i.e. do not de- 
pend on the azimuth (f>, Eqs. (5) decouple into a 
set of 2 (iV + 1) equations. 

dtT"" = RnV^ + RaV^ + RmV^ 

+ riWlT^ + dMr)^drirT"^), 
dtS"^ = RnU^ + RcU^ + RmU^ 

+ n^ls^- (23) 

where is the Laplacian operator (21), in which 
d"^ /d'^cj) is replaced by —rm?. Each of the m-modes 
can therefore be found from an independent set 
of equations. The functions C/™, 1^™,?/™, V^, 



W^, and V^, arc solutions of the m- transformed 
equations (14-19), such as: 

= [r • V X 1™ aB™ (r, 6») • e*'"'^ , (24) 

= [r-Vx Vx]™aB'"(r,^)-e'™^, (25) 

for the potentials C/™, V^". Analogous expressions 
can be found for the other potentials in Eqs. (14- 
19). The operators [r • Vx]", [r • V x Vx]™ and 
[— L^]™ are the m-transformed operators, where 
the derivatives with respect to ^ become multipli- 
cations by the complex factor im. 

We solve the set of equations (23) numerically. 
At each time step, we solve Eqs. (24,25) and the 
corresponding equations for the u x B term for the 
auxiliary potentials ,u]S,u^ . Equa- 

tions (23) are then advanced in time with a third- 
order Rungc Kutta method. Spatial derivatives 
are second-order accurate in a regular [r, 9] mesh. 
The models below have the resolution of 80x160 
mesh points. 

2.3. Coupling of the Modes 

When the assumption of axisymmetry for a is 
relaxed, the set of dynamo equations (23), for 
each m, are no longer independent. But it is 
still possible to use Eqs. (24,25), where the product 
a{r, ^)B'" is replaced by the convolution 

oo 

(aB)™= J2 a^{r,9)B'^-^{r,9). (26) 

j=-oo 

The coefficients with negative indices are equiv- 
alent to the complex conjugate of the coefficient 
with the positive indices. 

We investigate in detail the simplest situation 
of coupling of m = and m =1 modes, retaining 
only the lowest m = and m = 1 terms in the 
expansion (26). The contribution to the m = 1 
mode thus becomes: 

(aB)^ =a°B^ +ea^B°. (27) 

The non-axisymmetric part of a (a^) thus intro- 
duces a linear coupling between the m = 1 and 
the m = modes. The amplitude of the resulting 
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non-axisymmetric field depends on the value of e. 
In the results displayed below, e = 1 is used in 
model Ml and e = 0.2 in models M2 and M3. For 
an aQ dynamo, when ^ Ra, Eq. (27) means 
that the poloidal m = 1 component of the field 
is driven by the toroidal m = component. The 
toroidal m = 0, m ~ 1 modes are related respec- 
tively to the poloidal m = 0, m = 1 modes by the 
differential rotation. 

2.4. Boundary conditions 

The conditions on the rotation axis 9 = are 
different for different values of m. For even m, 
any function maps into itself when (f> ^ (p + t: 
and therefore on the axis the derivative in the 9- 
direction perpendicular to the axis must vanish. 
For odd m every function subjected to this map- 
ping is multiplied by —1 and therefore must vanish 
on the axis {9 = 0) to satisfy continuity. Hence, 



95 _ ST _ 

'de~'d9~' 

S = T = 0, 



m = 0, 



m = 1. 



(28) 
(29) 



We consider a spherical shell ri„ < Rq < 1 with 
internal radius ri„ = 0.4i?Q, where we define our 
potentials S, T to be zero. The physical boundary 
between the core and the convection zone is identi- 
fied by the sharp change in electrical resistivity at 
0.69i?o. The computational domain extends fur- 
ther down, to O.4i?0. We allow enough time for 
the diffusion and differential rotation to expel the 
flux from the low-conductivity region. 

On the outer boundary, r = R, we adopt the 
potential (vacuum) field conditions: 



0, 



dsi 

dr 



0, 



(30) 



where 5™ are the coefficients of the decomposition 
into spherical harmonics of degree £ and order m 
of the scalar potential S at the surface. For the 
numerical scheme, which is second order in the r- 
direction, we only need to solve this equation at 
two r-layers which are next to the boundary. 

2.5. Meridional Circulation 

Meridional circulation has long been observed 
on the surface of the Sun with a variety of tech- 
niques from sunspot tracking (Tuominen 1941) to 



Dopplcr and helioseiseismic measures (LaBonte & Howard 1982; 
Hathaway et al. 1996; Haber et al. 2002). 

The patterns correspond to a single-cell large- 
scale flow proceeding from the equator to the pole 
with an amplitude at the surface of about 20 m/s. 
Recent measurements indicate, however, that the 
flow might in fact have a more complex nature, 
even reversing at times and showing a non-trivial 
vertical structure (Haber et al. 2002). 

With current techniques, the upper 15% of the 
solar surface can be probed (Braun & Fan 1998). 
The return flow is not observed although some 
progress has been made in this direction (Giles et al. 1997). 
In particular, the depth into the convection zone 
to which it penetrates is not known. Hathaway 
et al. (2003) examine the drift of the centroid of 
the sunspot area toward the equator in each hemi- 
sphere from 1874 to 2002 and find that these ob- 
servations are consistent with a meridional coun- 
terflow deep within the Sun as the primary driver 
of the migration toward the equator and the pe- 
riod associated with the sunspot cycle. 

In the context of an axisymmetric dynamo, 
meridional circulation has been extensively stud- 
ied (Roberts & Stix 1972; Wang et al. 1991; Choudhuri et al. 1995; 
Durncy 1995; Dikpati & Oilman 2001; Nandy & Choudhuri 2002). 
For the magnetic fields to be effectively trans- 
ported by such a small amplitude flow, the elec- 
trical resistivity has to be low, a condition that is 
fulfilled below the convection zone, where resistiv- 
ity and turbulence dramatically decrease. Hence, 
one can expect that a deep flow, sinking below the 
convection zone, will be more effective in acting 
upon the field (Nandy & Choudhuri 2002). How- 
ever, it has not been established that such a flow 
can realistically exist. We will discuss two models, 
see Fig. 3. 

2.6. Initial Conditions 

We select an initial field of antisymmetric (dipo- 
lar) symmetry with respect to the equator by 
choosing the following potentials for m = and 

m = I modes: 

S{t = 0) = S°, S\t = 0) = Sl, 

T{t = 0) = T^, T\t = 0) = Tl (31) 

Here, 5™ indicates a function proportional to the 
spherical harmonic of order I and degree m. No 
symmetry with respect to the equator is imposed 



5 



on the time-dependent solution, and the whole do- 
main G [0, it] is therefore included. 

The velocity field u in the induction equation 
(1) is given. The induction equation is thus lin- 
ear in the magnetic field B and its solutions are 
exponentially growing or decaying in time. A sta- 
tionary state is often achieved by introducing a 
quenching of the a-efltect when the field intensities 
exceed some level, for a discussion see for instance 
(Blackman & Field 2002). We do not prescribe 
here any quenching mechanism because it would 
introduce a further coupling between modes, mak- 
ing the interpretation of the results more com- 
plex. It would also add an element of arbitrari- 
ness, which is not necessary for our purposes. We 
will therefore be looking at solutions that are close 
to marginal excitation. 

3. Results 

We consider first the dynamo without merid- 
ional circulation (um =0). In Table 1 the growth 
rate 7, the adimensional cycle period T and the ra- 
tio of the cycle period to the core rotation period 
T /Trot as explained in Section 2, are given. The 
three models considered are labelled as Ml, M2, 
M3. See also the Appendix for a more detailed 
description. 

3.1. Growth rates and cycle periods 

The time evolution of the toroidal field poten- 
tial for the Ml case is shown in Fig. 4. A steady 
exponential growth is reached after one diffusion 
time. We follow the; run for approximately five 
diffusion times. However, already after t ^ 1, the 
solution appears to be stable. As seen in Fig. 4, 
the solution is oscillating. Each cycle lasts 179 ro- 
tations of the body being modelled. The ratio of 
the cycle period to the rotation rate, {T /Trot) is 
independent of the units adopted. This ratio dis- 
criminates between solutions that have a solar-like 
period and those which do not. For model Ml, the 
ratio T/Trot is 179 and is close to the solar case of 
149. 

Model M3, not shown, behaves similarly to 
model Ml. It has a comparable oscillation period. 
The a-effcct in these two models have a broad ra- 
dial distribution within the convection zone which 
overlaps with the shear layer, see Fig. 2. 

In contrast, model's M2 behaviour is markedly 



different than the previous two, especially regard- 
ing to the oscillation period T, which is ten times 
smaller than the cases Ml and M3. T/Trot of 
20.2 shows that M2 does not represent the solar 
case correctly. M2 is characterised by the a-effect 
which has no overlap with the shear layer and this 
is probably the cause of its failure. 

3.2. Phase relations between the axisym- 
metric and non-axisymmetric modes 

In Fig. 5 we show the time evolution of the ra- 
dial and azimuthal components of the magnetic 
field during one cycle in model M3. The axisym- 
mctric dipolar-type field at the pole reverses its 
sign when the (total) azimuthal field close to the 
tachocline is maximal. In other words, there is a 
phase shift of 7r/2 between these two field com- 
ponents. We evaluate the toroidal fields at 30° 
of latitude, where the non-axisymmetric field is 
concentrated. For the plot of the to = 1 compo- 
nent, which has the cos<^ dependence, we select 
the meridional plane defined by = 0, where the 
field is maximum. A similar relation is observed 
on the Sun: the axisymmetric dipole reverses close 
to the solar maximum. Analogous calculations for 
models Ml and M2 yield 7r/6 and — 3/47r, respec- 
tively, see Table 2. We may therefore conclude 
that the phase relation between the field compo- 
nents is dependent on the radial distribution of 
the a-effect. This can help to select the models 
according to the phase relations between the field 
components, as originally indicated by Stix (1976). 

Another phase relation that can be observed 
is the one between the surface axisymmetric and 
non-axisymmetric poloidal fields is shown in Fig- 
ure 9 (Ruzmaikin et al. 2001). Note that, accord- 
ing to their definitions, the poloidal potentials 
have the same phase as the poloidal fields. For 
the purpose of determining the phase difference 
between the modes of the poloidal field, we can 
therefore use the potentials, see Figs. 6-8. 

The time dependence of the m = and m = 
1 poloidal potentials is shown in Figs. 6-8. The 
phase relations can be read from the plot and are 
fisted in Table 2. 
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3.3. Spatial distribution of the fields. Lo- 
calisation of the non-£Lxisymmetric 
modes 

The poloidal field is found to be predominantly 

antisymmetric (odd) with respect to the equator 
(dipolar-type) in all three models. However, this 
symmetry is not pure and a symmetric component 
is present as well, see Fig. 10,11. We recall that we 
do not impose an equatorial symmetry conditions 
on the solution and solve the equations for the 
whole domain 9 G [0,7r]. The north-south symme- 
try is broken because we couple the axisymmet- 
ric mode, which is odd across the equator (like 
polar dipole), with the non-axisymmetric mode, 
which is even across the equator (like equatorial 
dipole). This north-south asymmetry does not dis- 
agree with observations. Fig. 12 shows the time 
evolution of the signed azimuthal magnetic field 
for the m = 1 mode on the meridional plane (j) = Q 
during half of a cycle. In contrast to the axisym- 
metric mode, the non-axisymmetric field is highly 
concentrated in the lower part of the tachocline 
and, to lesser extent, near the surface. It is also 
localised around 30° latitude, see also Table (3). 

In order to show the robustness of the spa- 
tial intensity distribution of the fields throughout 
a cycle we display the distribution of the abso- 
lute values of the the axisymmetric (m =0) and 
non-axisymmetric (m = 1) toroidal fields inside 
the convection zone, integrated over a solar cycle 
(Figs. 13,14,15). In Table (3) we show the loca- 
tions of the maxima of the toroidal m = and 
m = 1 components, r\j, 6^, r\f and 9\f. 

The physical reason for these localisations can 
be understood. As it is well known from MHD 
theory (Mofffat 1978; Krause & Radler 1980), the 
differential rotation affects the axisymmetric and 
non-axisymmetric fields in different ways. Thus, 
an axisymmetric poloidal field in a differentially 
rotating region is converted into an axisymmet- 
ric toroidal field. But a non-axisymmetric field 
is rapidly (within a few rotations) excluded from 
the strongly differentially rotating high resistivity 
region. This happens because the rotation tends 
to twist the field lines perpendicular to the axis 
of rotation in such a way that any two adjacent 
lines would have an opposite direction. This leads 
to effective rcconncction of these field lines, ex- 
cluding the field from the differentially rotating 



region even in the case of a very small electrical 
resistivity. Hence, differential rotation is a mech- 
anism that tends to destroy any deviation from 
axial symmetry. In the solar convection zone, the 
field survives only in the regions where the dif- 
ferential rotation is weak (around 30°) and where 
the resistivity is small, i.e. in the lower part of the 
tachocline. 

The localisation of non-axisymmetric magnetic 
field suggests a simple explanation of why ac- 
tive regions emerging on the solar surface have a 
tendency to cluster on preferred longitudes. Ac- 
tive regions, such as sunspots, are the result of 
buoyancy of bipolar magnetic loops (flux ropes) 
through the solar surface, (Parker 1979). The 
flux ropes formed at the base of the convection 
zone emerge when their magnetic fleld exceeds 
the threshold for buoyancy. A non-axisymmetric 
enhancement of the underlying magnetic field at 
that location results in the clustering of active re- 
gions as shown in a heuristic model (Ruzmaikin 
1998,2001). The m = 1 mode of the toroidal 
magnetic field superimposed on the axisymmet- 
ric mode produces a localised field maximum 
( "hump" ) near the maximum of sin (j). The humps, 
produced by non-axisymmetric fields, are unsta- 
ble when their field strength reaches 4 x 10^ G 
(Caligari et al. 1995; 1998). Because the growth 
of the underlying non-axisymmetric mode is af- 
fected by the continuous stretching by differential 
rotation, the mode is after some time destroyed. 
The stretching time scale of several months is in 
agreement with the observed life-time of the clus- 
ters (Gaizauskas et al. 1983). 

3.4. Rotation Rate and Preferred Longi- 
tudes 

As we mentioned in the Introduction, a per- 
sistent, periodic magnetic pattern in the so- 
lar wind has been observed, with a fixed ro- 
tation rate (Neugebauer et al. 2000). A pe- 
riod consistent with this rate was also found 
from the surface fields (Ruzmaikin et al. 2001; 
Henney & Harvey 2003). If preferred longitudes 
are associated with the non-axisymmetric modes 
of the dynamo, then the observed pattern would 
originate in the poloidal non-axisymmetric com- 
ponent of the field. 

The longitudinal distribution of the non-axisymmetric 
modes of the dynamo has the form of a wave 
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propagating along the azimuth (f> with a rate de- 
termined by the solution of the dynamo equations. 
We determine that rate, , by following the time 
evolution of the surface poloidal (radial) m = 1 
mode, Fig. 16. The results are given in Table 1 for 
the three models considered. We represent these 
rates as horizontal dashed lines, along with the the 
radial distribution of the rotation rate, in Fig. 1. 

The rotation rate of the non-axisymmetric m = 
1 mode is comparable to that of the surface rota- 
tion at 30° latitude for the Ml and M3 models. 
This is also the rotation rate of the core. For the 
M2 model, a slower rate is found. Recent Ulysses 
observations of the surface magnetic fields during 
Ulysses' fast scan around solar maximum, show 
a rotation rate of the non-axisymmetric dipole 
of the Sun close to the core rotation rate, see 
(Jones et al. 2003). 

3.5. Meridional Circulation 

The main goal of this paper is to study the 
excitation of the non-axisymmetric mode un- 
der otherwise simple conditions. For this pur- 
pose we did not include the meridional circu- 
lation in the models described above. Merid- 
ional circulation may play an important role 
in the overall dynamo, as indicated by other 
studies (Roberts & Stix 1972; Wang et al. 1991; 
Choudhuri et al. 1995; Durney 1995; Dikpati & Oilman 
Nandy & Choudhuri 2002). Here we report only 
on a preliminary analysis of two simple cases of 
meridional circulation for the Ml distribution of 
a. The cases are (see Appendix): (1) shallow and 
(2) deep meridional circulation, see Fig. 3. In 
order to make the meridional transport more ef- 
fective compared to diffusion in these calculations 
we decrease the ratio of the core resistivity to the 
resistivity in the convection zone to 1/200. 

The resulting field distributions two cases are 
displayed in Figs. 17,18. Comparing these distri- 
butions with the case without meridional circula- 
tion. Fig. 13, we see that the field tends to be de- 
pleted in the upper part of the convection zone and 
intensified and stretched in the tachocline. This is 
especially true for the case with deep meridional 
circulation. 

The field is pushed out of the convection zone 

into the tachocline, apparently under the action of 
the radial component of the circulation. The flow. 



however, does not seem to be effective enough in 
the latitudinal direction to significantly change the 
toroidal field distribution, although there is a shift 
of its maximum below the 30° line (see Table 4). 

A more significant change occurs in the cycle 
period of the solutions, which nearly doubles with 
an amplitude of the flow |um = 5| in comparison 
with the case |um| = 0. 

Moreover, with |um > 5|, the symmetry of the 
solution changes from dipolar (antisymmetric with 
respect to the Equator) , to quadrupolar (symmet- 
ric). For um = 10 a steady non-oscillatory solu- 
tion takes over. 

Discussion 

Solar observations have shown the importance 
of the non-axisymmetric component of the mean 
magnetic fields in the formation and evolution of 
solar activity. Solar dynamo models have so far 
failed to address the non-axisymmetric features 
of the solar cycle, most notably the clustering of 
magnetic activity and preferred longitudes. 

We have developed a new kinematic non- 
axisymmetric mean-field dynamo model of the Sun 
in spherical geometry, which incorporates the solar 
rotation as reconstructed from helioseismic data 
and a model for the meridional circulation. We 
2(^f.ompose the magnetic field into its toroidal and 
poloidal components. This geometry allows waves 
to propagate azimuthally at different rates for dif- 
ferent modes, contrast to a periodic box where 
only one rate can be allowed. We determine the 
rotation rate of the first non-axisymmetric mode. 

We have considered for simplicity only the low- 
est azimuthal modes, m = and m = 1. By the 
action of a non-axisymmetric a-effect, the modes 
are coupled together, and the m = I mode shares 
the same cyclic behaviour as the axisymmetric 
mode m = 0. 

We have examined three different radial distri- 
butions of the mean helicity (a-effect), two with a 
operating near the tachocline and one only near 
the surface. When the meridional flow is ne- 
glected, we found that all produced a cycle with al- 
ternating m = and m = 1 modes. The two mod- 
els with a near the tachocline both give the ratio 
T/Trot of the cycle duration to the rotation rate 
within 20% of the observed value of 147, whereas 
the case where a was concentrated near the sur- 
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face failed to give agreement indicating that this 
was not a valid solar model. 

The radial profile of a strongly influences the 
phase relations between the alternating axisym- 
metric dipole and the toroidal field intensity at 
the base of the convection zone. 

We determined the rotation rate of the near- 
surface non-axisymmetric radial field in the three 
models. The two models with a close to the 
tachocline had the same rate of rotation, slightly 
faster than the core rotation rate, in good agree- 
ment with observations. This rate was found to be 
the same through several cycles. Again, the model 
with near-surface a was different, with a rotation 
rate slower than the core rotation. 

The spatial distribution of the non-axisymmetric 
toroidal fields was found to be concentrated 
around 30° latitude, where the radial gradient 
of the differential rotation vanishes. [The lati- 
tude, at which the axisymmctric mode is con- 
centrated is higher. However, in axisymmct- 
ric dynamos, it is possible to vary the merid- 
ional circulation and a latitudinal distribution 
of a to lower the latitudinal distribution of 
toroidal fields (Ruediger & Brandenburg 1995; 
Dikpati & Oilman 2001). ] When a was finite 
near the tachocline, the non-axisymmetric toroidal 
field had a strong component there. A near surface 
concentration was also found for all three models. 

A preliminary study of meridional circulation 
indicated that the overall field distribution is not 
strongly affected, particularly for the toroidal non- 
axisymmetric component of the field. A bigger 
effect was found in the cycle period, which be- 
comes longer at higher flow intensities, and the 
equatorial symmetry of the solutions, which be- 
comes quadrupolar (symmetric) when flow inten- 
sity is very high. 

In a forthcoming study we will investigate 
the role of meridional circulation in the non- 
axisymmetric dynamo in more details . 
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A. Description of sources 
A.l. Rotation 

The differential rotation is approximated by 

^(r, ^) = ^ (^^30 (1 - /(r)) + nsfir)^ , /(r) = ^ (^1 + tanh (^^) ) • (Al) 

The parameters rn = 0.69 and 6n = 0.05 define the location and the width of the tachocline. Here, Cls is 
the surface rotation profile, 

Os (6') = 1 + ai cos^ e + a2 cos^ 0, ( A2) 

with ai = -0.126, 02 = -0.159, , and 6 is the co-latitude (Moss 1999; Charbonneau et al. 1999). The 
equatorial rotation rate is thus given by 

fie, = ^ (A3) 
"30 

where O30 = Os(30°). In all runs we use and Oo = 4.9 x 10*, in non-dimensional units. 

A.2. DifFusivity 

The turbulent diffusivity profile is given by 

77(r)=77c(l-/(r))+%/(r), /(r) = ^1 +tanh (^^^^^ . (A4) 

Here, s m = 0.69, (5^ = 0.04; r]c and rjo are the diffusivity values in the core and in the convection zone 
respectively. We use ?7o/??c = 50 and r]o/r]c = 200. 

A. 3. a-effect 

The radial dependences of the a-effect for the three cases considered, are: 



»{r) ="o^y- + tanh (A5) 

with ao = 0.9, = 0.7, 6a = 0.03, in model Ml, and ao = 24.2, Va = 0.9, Sa = 0.03 in model M2; 

air)=ao - M +tanh^-^ j ■- f 1-tanh ^ j (A6) 

with ao = 1.9, = 0.75, Sa = 0.04, for model M3. The values of ao have been chosen so that the dynamo 
is supercritical but close to marginal. 

A. 4. Meridional circulation 

We consider a simple 2-cell flow where the surface flow amplitude is of the order of that observed. Density 
stratification in the convection zone is very high, therefore mass continuity has to be taken into account. 
We specify a simple fit to reproduce the density from a solar model (Christensen-Dalsgaard et al. 1996), 
between 0.60i?© and 0.95i?o. 

p = ae-'^'' (A7) 
where a = 247 and /3 = 10.1. We impose mass continuity 

Vpu = (A8) 
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by introducing a poloidal potential ^ for the mass flow: 

pu = V X V X f*(r, 6). (A9) 

The analogous Stokes stream function is — rsin^ dg^. Streamlines for this flow are contours where such a 
function is constant. Fig. 3 shows the two cases of the meridional circulation where the flow does or does 
not penetrate inside the shear region (tachocline). Those profiles result from a simple power-law expression 
(Roberts & Stix 1972) for the radial part of the potential: 

2 

(1 - rbj" 

^{r) = r <n, (AlO) 

where = 0.45 for the deep flow case and rb = 0.65 for the shallower case. We specify the latitudinal 
dependence to be proportional to 

P2,o(cos0) - Scos^e - 1 (All) 

where Pi^m is the associated Legendre polynomial of order I = 2 and degree m = 0. This splits the domain 
into two counter-rotating cells, each poleward of the equator. 
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Fig. 1. — Distribution of the internal solar rota- 
tion used in our dynamo modeling. Dashed curves 
show the rotation rates of the non-axisymmetric 
mode, see Table (1). The upper curve at 442nHz, 
just above the core rotation rate, is for models Ml 
and M3. The lower curve, at 433nHz, is for model 
M2. 
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Fig. 3. — Meridional circulation. Flow lines for 
the two cases considered. Solid lines indicate anti- 
clockwise flow and dotted lines clockwise. The 
dashed line at 0.69/?© indicates the centre of the 
tachocline. 
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Fig. 2. — Radial distributions of the sources of 
magnetic fleld generation: a, for three different 
models, the radial gradient of the rotation rate at 
the equator and turbulent diffusivity rj. Values are 
not to scale; for actual values see Appendix. 
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Fig. 4. — Long-term time evolution of the ampli- 
tude of the toroidal field potentials for model Ml. 
Top panel is Tq and lower panel is Tq. The diffu- 
sion time is taken as the time unit. 
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Fig. 5. — Phase relation between magnetic field 
components in model M3. The time evolution dur- 
ing a cycle of magnetic field components is plotted. 
The solid line shows the intensity of the azimuthal 
component of the total magnetic field, B° + B^, 
in the meridional plane (p = 0, calculated at the 
tachocline, r = O.69i?0, at 30° latitude. The dot- 
dashed line shows the intensity of the surface ra- 
dial field component B° close to the north pole, at 
80° latitude. The phase shift between the fields is 
7r/2. Amplitudes of both toroidal and radial fields 
are normalised to unity for convenience. 




T, 
T, 



4.929 4.939 4.949 4.959 4.969 

t (R'/t,) 

Fig. 6. Phase relations between poloidal and 
toroidal components of the axisymmetric and non- 
axisymmetric magnetic fields for model Ml. The 
integrated square moduli of the scalar potentials 
S,T, over the whole domain are plotted versus 
time. Values are not to scale. Tq and are 
the toroidal and poloidal potentials for the m = 
mode, Ti and Pi are the toroidal and poloidal po- 
tentials for the m = 1 mode. The Si mode is cou- 
pled to the To mode (sec Eq. 27). The observed 
phase relations for the Sun's field components are 
shown in Fig. 9 for comparison. 
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Fig. 7. — Same as in Fig. 6 for model M2. 
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Fig. 8. — Same as in Fig. 6 for model M3. 
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Fig. 9. — Observed phase relations between 
poloidal components of the axisymmetric (dashed 
Hne) and non-axisymmetric (solid line) magnetic 
field (Ruzmaikin et al. 2001). 
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Fig. 10. — Distribution of the axisymmetric (rn = 0) poloidal field over half a cycle. Time proceeds from 

left to right and from top to bottom and panels arc equally spaced in time. Poloidal field lines are shown. 
Solid lines have a counter-clock wise direction. Dotted lines mean opposite direction (model Ml). 
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Fig. 11. — Butterfly diagram for the rn = (dot-dashed contours) and m = 1 (solid contours) modes of 
the toroidal magnetic field at two different depths: above the shear layer, O.74i?0 (a), and in the convection 
zone, at O.85i?0 (b). Model Ml is displayed here. The diagram shows the latitudinal distribution of the 
modulus of the toroidal magnetic field as a function of time, during five half-cycles. The time is measured 
in unit of the diffusion time. 
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Fig. 12. — Distribution of the non-axisymmetric {jiiy^ 1) toroidal field through half a cycle, until the field 
has changed its polarity. Model Ml. The contour plot show the intensity of the azimuthal component of 
the non-axisymmetric field B^. The dashed line marks the centre of the tachocline, at r = O.69i?0, which is 
also the location of the turbulent resistivity drop. 




Fig. 13. — Distribution of the axisymmetric (m = 0) and non-axisymmetric (m = 1) toroidal fields inside 
the convection zone in model Ml. The modulus of the azimuthal component of the magnetic field is shown. 
Field distribution is integrated over time during a cycle. Contour lines concentrate into region where the 
field is more intense and resides for longer times. For the m = 1 component, the modulus squared of the 
field is averaged over the azimuthal direction The dashed line marks the centre of the tachocline, at 
r = O.69i?0 which is also the location of the turbulent resistivity drop. 




Fig. 14. — Same as in Fig. 13 for model M2. 
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Fig. 16. — Left panel: Near-surface distribution of the radial component of the non-axisynimetric magnetic 
field at fixed time. White indicates outwardly directed field, black inwardly. Right panel: amplitude of the 
radial field versus longitude at the latitudinal cut indicated by the line in left panel. Two curves are shown, 
obtained at different times during a cycle. The phase diflFcrencc between the maxima divided by the time 
lag between the two profiles returns the rotation rate of the m = 1 mode. 
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Fig. 18. — Same as in Fig. 13 for model with deep meridional circulation. The flow lines are superimposed 
as dot-dashed lines. 
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Model Ml M2 M3 



7 0.35 4.15 1.4 

T 0.022 0.0026 0.018 

T/Trot 179 20.2 137 

9} (nHz) 442 433 442 



Table 1: Time evolution. 



Model 




Ml 


M2 


M3 


A<^ {B^ - 




l/67r 


-3/47r 


l/27r 


{So - 




l.lTT 


1.47r 


l.lTT 



Table 2: Phase relations. 



Model Ml M2 M3 

070 095 0.72 ; 0.90 

6'^(lat) 57 50 60 ; 52 

rlj 0.67 ; 0.93 0.90 ; 0.97 0.78 ; 0.93 

6lj^(lat) 30 ; 33 42 ; 41 32 ; 32 



Table 3: Location of magnetic field maxima. 



Model MIS MID 

't^ 080 065 

6'^(lat) 55 59 

r\t 0.67 ; 0.90 0.65 ; 0.87 

6>j^^(lat) 31 ; 34 30 ; 35 



Table 4: Location of magnetic field maxima. 
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